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ABSTRACT 


A  ducted  rotor  was  tested  in  the  36-inch  water  tunnel  at  David  Taylor 
Model  Basin,  Naval  Surface  Warfare  Center/Carderock  Division.  At  the  blade 
tip,  both  leakage  and  trailing  edge  vortices  were  observed  through  cavitation 
inception.  Detailed flow  patterns  in  the  vortices  have  been  mapped  with  LDV 
(LaserDoppler  Velocity)  technique.  Numerical  calculations,  based  on  a 
Reynolds-averaged  Navier-Stokes  formulation,  were  performed  to  simulate  the 
observed  vortices.  To  reduce  the  numerical  diffusion,  a  fifth-order-accurate 
upwind  differencing  scheme  was  implemented  to  evaluate  the  flux  difference 
across  the  cell  surface.  The  results  were  compared  with  the  experimental  data. 

ADMINISTRATIVE  INFORMATION 

This  work  was  sponsored  by  the  Advanced  Submarine  R&D  office  (93R),  Naval  Sea 
Systems  Command  (NAVSEA)  and  administrated  by  the  Advanced  Propul  sor  Development 
Program  Office  (Code  508),  Naval  Surface  Warfare  Center,  Carderock  Division  (NSWCCD).  The 
work  was  conducted  by  the  Hydromechanics  Directorate,  Propulsion  and  Fluid  Systems 
Department  (Code  5400),  NSWCCD,  under  work  unit  number  1-5080-262-52. 

INTRODUCTION 

This  report  is  about  a  numerical  simulation  of  tip  vortices  of  a  rotor  operating  in  a  duct. 
The  results  of  the  simulation  are  to  be  compared  with  the  test  data  of  rotor  P-5206.  The  rotor  has 
three  blades,  and  is  designed  to  operate  in  the  36-inch  (91 .44  cm)  water  tunnel  for  the 
investigation  of  the  tip  vortices  cavitation  of  a  ducted-rotor.  The  water  tunnel  is  of  an  open  jet 
type,  the  flow  is  discharged  into  a  90.0  inch  (228.60  cm)  plenum  chamber  through  the  nozzle. 

The  pressure  in  the  tunnel  can  be  adjusted  so  that  the  cavitation  inception  can  be  investigated. 

The  rotor  is  positioned  such  that  its  blade  tip  trailing  edge  is  upstream  of  the  nozzle  of  the  open 
jet.  Figure  1  shows  the  positioning  of  the  rotor  in  the  36-inch  water  tunnel.  The  cylindrical  tunnel 
wall  serves  as  the  duct  for  the  rotor.  The  detailed  rotor  geometry  and  the  test  conditions  are  given 
in  Reference  [1],  only  some  relevant  geometrical  dimensions  will  be  mentioned  briefly  in  the 
following.  The  rotor  blade  has  a  constant  chord  length  of  15.0079  inches  (38.12  cm)  and 
maximum  section  thickness  of  1.50  inches  (3.81  cm).  The  diameters  of  the  rotor,  its  hub  and  the 
inner  surface  of  the  duct  are  33.475  inches  (85.03  cm),  13.846  inches  (35.17  cm)  and  34.0  inches 
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(86.36  cm),  respectively.  The  skew  angles  at  the  root  and  at  the  tip  are  -4  degree  and  30  degree, 
respectively.  The  ratio  of  the  tip  gap  to  the  maximum  section  thickness  is  0.175.  The  axial 
distance  between  the  trailing  edges  of  the  nozzle  and  the  blade  tip  is  6.27  inches  (15.93  cm).  The 
rotor  was  driven  from  upstream.  The  velocity  measurements  behind  the  rotor  were  conducted 
with  a  three-component  LDV  system.  The  system  consisted  of  two  TSI  model  9832  fiber  optic 
probes  rigidly  mounted  together  on  a  traverse  that  could  translate  in  the  axial,  vertical,  and 
horizontal  directions.  The  detailed  techniques  used  to  measure  and  analyze  the  flow  field  were 
presented  in  Reference  [2], 

Recently  two  numerical  calculations  about  the  tip  vortices  of  propeller  P-5206  were 
reported.  One  of  the  calculations,  based  on  a  finite-volume  approach  with  artificial 
compressibility  algorithm  [3],  was  by  Brewer  [4];  and  the  other,  based  on  a  finite-differencing 
approach  with  pressure  implicit  split  operator  algorithm,  was  by  Kim  [5].  The  former  calculation 
was  supported  by  an  unstructured  topology  while  the  latter  calculation  was  supported  by  a 
structured  topology  with  overset  capability  derived  from  CHIMERA  scheme  [6].  Both 
calculations  used  upwind  method  to  evaluate  the  convective  terms.  However,  in  the  former  case, 
the  differential  system  is  hyperbolic  and  the  upwind  direction  was  determined  by  the  signs  of  the 
eignvalue  of  the  system;  while  in  the  latter  case,  the  upwind  direction  was  determined  by  the  sign 
of  the  convecting  vectors,  since  the  differential  system  is  elliptic  in  nature.  Despite  of  the  distinct 
differences  in  algorithm  and  scheme,  both  calculations  yield  similar  and  excellent  results  in 
certain  aspects.  This  report  presents  the  results  of  a  third  calculation.  The  solution  scheme  is 
similar  to  Brewer’s  [4],  except  that  the  calculations  were  supported  with  a  structured  rather  than 
an  unstructured  topology.  Although  the  unstructured  topology  is  more  robust  in  generating  the 
supporting  grid  system  and  the  grid  can  be  easily  refined  at  where  it  is  needed  the  most,  it  has 
been  reported  that  the  flow  field  near  the  solid  boundary,  where  most  of  the  vorticity  is  generated, 
are  predicted  better  with  a  structured  grid  topology.  Perhaps,  it  is  because  the  structured  grid  can 
better  resolve  the  laminar  sub-layer  near  the  solid  surface.  Overset  topology  [5]  eliminates  many 
difficulties  associate  with  the  grid  generation  and  improves  grid  quality  about  complicated 
geometry;  however  interpolations  are  needed  in  order  to  communicate  between  the  various 
overset  grid  systems,  and  interpolations  often  introduce  errors. 

Present  calculations  are  complementary  with  the  two  mentioned  above.  The  solver  used 
is  the  UNCLE  code  [7],  and  the  grid  topology  is  structured  in  nature.  The  structured 
computational  grid  topology  used  here  is  designed  such  that  it  is  able  to  resolve  the  separated 
flow  near  the  blade  tip  and  it  can  easily  be  refined  to  support  the  passage  of  the  vertical  flow. 
Calculations  based  on  two  different  grid  sizes  were  performed.  The  first  calculation  was 
exploratory;  its  goal  is  to  obtain  an  assessment  about  the  fidelity  of  the  solutions,  such  as  the 
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pressure  and  the  helicity  (dot  product  of  vorticity  and  velocity  vectors)  distributions  about  the 
core  of  the  vortices.  Thrust  and  torque  of  the  rotor  blade  should  be  predicted  reasonably  well.  The 
goal  of  the  second  calculation  is  to  study  the  sensitivity  of  the  solutions  with  respect  to  the  grid 
refinement.  Based  on  the  trajectory  of  the  tip  vortices  of  the  first  calculation,  the  grid  was  then 
refined  to  produce  the  second  set  of  the  grid  for  the  refined  calculation.  The  grid  sizes  are  about 
3.0  and  5.5  millions  nodes  for  the  first  and  the  second  set,  respectively.  In  addition  to  the  standard 
MUSCL  scheme  (monotone  upstream-centered  scheme  for  conservation  laws)  provided  by 
UNCLE  code  [7],  a  third-order  and  a  fifth-order  weighted  MUSCL  schemes  are  implemented. 
This  is  an  effort  to  reduce  the  numerical  diffusion  in  the  vortex  core  region. 

The  calculated  flow  quantities  such  as  velocity  components,  helicity  are  compared  with 
the  measured  values  and  the  pressure  distribution  in  vortex  core  is  compared  with  previous 
computations. 


DESCRIPTION  OF  THE  EXPERIMENT 


The  detailed  description  and  setup  of  the  experimental  facility  and  the  offset  of  the  blade 
section  can  be  found  in  references  [1]  and  [2].  The  rotating  speed  of  the  rotor  was  500  rpm;  the 
tunnel  speed  was  22.85  f/s  (6.96  m/s).  The  Reynolds  number  based  on  tunnel  speed  and  the  rotor 
radius  was  3.019  million.  Flow  into  the  rotor  was  measured  in  the  duct  at  a  distance  of  0.361  rotor 
radius  upstream  of  the  rotor  stack  line.  The  inflow  was  almost  free  of  the  tangential  component, 
and  the  axial  component  of  the  inflow  is  shown  in  Figure  2.  The  thickness  of  the  boundary  layer 
near  the  duct  surface  is  much  larger  then  the  tip  gap  of  0.0 15R.  The  rotor  tip  is  well  imbedded  in 
the  inflow  boundary  layer.  Measurements  of  all  three  components  of  the  velocity  were  made  in 
the  tip-vortex  region  just  downstream  of  the  rotor  tip  at  103  different  x-r  planes.  Figure  3 
indicates  the  location  of  the  measuring  planes,  with  only  every  tenth  plane  is  shown.  The 
measuring  planes  are  normal  to  the  pitch  line  at  rotor  tip.  The  last  measuring  plane  is  about  half 
of  a  chord  length  from  the  rotor  trailing  edge.  Each  measuring  plane  was  divided  into  a  39  x  33 
mesh,  in  x  and  r  direction,  respectively.  At  each  mesh  point  approximately  100  measurements  of 
particle  velocity  were  collected. 

The  uncertainty  of  the  measurements  depends  on  the  sample  size,  the  flow  gradient  and 
the  turbulence  level.  It  also  depends  on  the  ability  to  control  the  tunnel  and  the  rotor  speed.  The 
maximum  uncertainty  of  the  axial  velocity  in  rotating  frame  occurs  in  the  vortex  core  and  in  the 
blade  wake,  they  are  10%  and  5%  of  the  tunnel  inflow  velocity  Vref ,  respectively.  The  uncertainty 
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in  the  ffeestream  is  less  than  1%.  The  uncertainty  in  the  inception  cavitation  number,  Oj  was 
estimated  at  10%. 


NUMERICAL  SCHEME 

The  numerical  solver  used  in  present  simulations  is  based  on  the  CFD  code  UNCLE 
developed  at  MSU.  It  is  based  on  a  three-dimensional  incompressible  Navier-Stokes  equation 
expressed  in  the  conservative  law  form.  The  equations  are  coupled  with  the  continuity  equation 
through  the  artificial  compressibility.  The  diffusive  terms  are  discretized  with  a  central 
differencing  scheme,  and  the  convective  terms  are  discretized  with  the  Roe’s  flux-difference 
splitting  scheme  [8].  For  the  steady  state  simulation,  a  first  order  Euler  backward  approximation 
is  used  to  represent  the  time  derivative.  To  accommodate  the  rotational  motion,  the  equations  of 
motion  are  first  formulated  in  a  rotating  coordinate  system  attached  to  the  moving  body,  but  are 
then  recast  in  the  absolute-flow  variables.  Such  formulation  is  popular  in  turbomachinery 
applications  [9,10].  It  is  believed  that  such  formulation  allows  more  accurate  calculation  of  flux 
in  the  finite  volume  method  [9]. 

The  Navier-Stokes  equations  of  the  three-dimensional,  time  dependent,  incompressible 
flow  can  be  written  in  a  curvilinear  coordinate  and  in  a  strongly  conservative  law  form  as 


+  (E-Ev)e+(F-Fv)11+(G-Gv)?=0 


where  Q  =  (p,  u,  v,  w) T 


sx  1l  fk 

‘E' 

J  J  J 

V 

and 

F 

__ 

Ik  Ik 

F* 

J  J  J 

* 

_F_ 

k  k  k 

G 

J  J  J 

with  the  inviscid  fluxes 

E*  =  (pu,u2  +p,uv,uw)T 
F*  =  (pv,uv,v2  +p,vw)T 


(1) 


(2) 


(3) 
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w,uw,vw,w2  + 


and  shear  fluxes 


EV,FV,G 


V 


P  is  the  artificial  compressibility  which  connects  the  continuity  equation  with  the  momentum 

equation,  J  is  the  transformation  Jacobian  which  relates  the  physical  domain  to  the  computational 
domain,  and  £x ,  £y  , £z ,  t|x ,  r|y ,  r|z ,  £x ,  ^y ,  £z  are  the  transformation  metrics. 


In  order  to  apply  Roe’s  concept  of  flux  differencing  to  the  convective  and  pressure 

differentiation  terms,  the  Reynolds  number  is  assumed  high  and  approaching  infinite.  The  shear 
fluxes  (Ev ,  Fv ,  Gv )  are  then  relatively  small  in  comparison  with  their  corresponding  inviscid 

fluxes  (E,  F,  G).  It  can  then  be  assumed  that  the  waves  that  travel  across  a  cell  face  are  normal  to 


the  face.  Subsequently  the  three-dimensional  system  of  partial  differential  equations  is  spatially 
split  into  three  quasi  one-dimensional  conservative  laws 


+  He  =0, 


(4) 


W  Vt 


where  0  is  £,  q.  or  Q,  and  H  is  E,  F,  or  G.  Defining  computational  cells  with  their  centroids  at 

m  (m=i,  j,  or  k),  and  their  cell  interfaces  at  m±l/2,  an  implicit  discrete  approximation  to  Equation 
4  can  be  written  as 


70 

+  DA 

m  m 

AQ"  +  AmHn 


=  0, 


(5) 


where  x  is  time  step  size,  AQn  =  Qn+1  -  Q"  and  AmH  =  [Hm+1/2  —  Hm_1/2]/A  0  and 
D  =  3H/3Q  the  flux  Jacobian.  For  the  construction  of  an  approximate  Riemann  solver  for  the 
initial  value  problem  in  Equation  5,  the  flux  at  the  cell  interface  m+1/2  that  separates  the  right 
from  the  left  traveling  wave,  can  be  define  as 


Hm+i/2  =H„  +(AH”+1/2)n, 


(6) 


or 


H 


n 

m+1/2 


=  Hnm+1-(AH;+I/2)\ 


(7) 
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or  the  simple  average  of  Equations  6  and  7.  At  present  calculation,  Equation  6  is  used,  since  it  has 
less  operation  count.  In  Equations  6  and  7,  the  superscript  -  and  +  of  the  flux  differencing  term 
represents  the  right  and  the  left  traveling  wave  across  the  interface  m+'A,  respectively.  That  is 

(AH;„J=[RA-R-'(Ql,-Q^)f.  (8) 

and 

(AH;.,n)‘=[RA*R-,(QLl-Q|;)f,  (9) 

with  R,  R'1  represent  right  and  left  eignvector  of  the  eigen  system  of  flux  Jacobian  D, 
respectively,  A+ ,  A-  represent  positive  and  negative  eigen  value  of  the  system,  respectively,  and 
Qn,+i  >  Qm  represent  the  dependent  variable  to  the  right  and  the  left  of  the  interface,  respectively. 
For  the  first  order  approximation,  Q^+]  =  Qm+1  and  QTni  =  Qm  . 

To  evaluate  the  Equation  6  numerically,  a  Roe’s  [8]  mean  value  matrix 
Dm+l/2  =  D(Qm,Qm+1)  is  introduced  with  the  following  properties; 

and 

D(Qm,Qm+1)  =  Dm,  if  Qm  — >  Qmn  •  (11) 

The  conditions  in  Equations  10  and  1 1  are  to  ensure  that  the  solution  to  Equation  5  is 
conservative  in  space.  For  the  incompressible  flow,  the  only  state  of  the  dependent  variable  Q  that 
satisfies  these  properties  is  the  simple  average  of  theQ^+1 ,  andQ^  .  Accordingly,  the  first  term 
on  the  right  hand  side  of  the  Equation  6,  H”  ,  is  computed  with  the  metrics  evaluated  at  the 
interface  m+'A  and  the  dependent  variable  Q  evaluated  at  centroid  m,  that  is  Q=Q^  •  For  the 
second  term  on  the  right  hand  side  of  the  Equation  6,  (A  Hm+1/2  )n ,  the  eignvectors  and  the 
eignvalues  R,  R  1 ,  A  ,  are  evaluated  with  the  metrics  and  the  dependent  variable  defined  at  the 
interface  m+'A,  that  is  Q=  \  (Q*+1  +  Qm  )  • 

To  enhance  the  accuracy  of  the  solution,  Van  Leer’s  [11]  MUSCL  (monotone  upstream- 
centered  scheme  for  conservation  laws)  approach  is  adopted.  The  dependent  variables  to  the  right 
and  to  the  left  of  the  interface  of  the  cell  are  reconstructed  as 
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QL,„ -Qm -i  [(i-kXQ„.2  -Q_.)-(i  +  kXQ«.  -q.)]. 


(12) 


and 

QL,«  =  Qm  +  f[(»  -  kXQ.  -  Q„-,)+  (i + kXQ„„  -  Q  J  •  (13) 

If  (p  =0,  the  Equations  12  and  13  recover  the  first  order  approximation.  If  q>  =1  and  K=  -1,  the 
approximation  is  second  order.  If  (p  =1  and  k=1/3,  the  approximation  is  third  order.  An 

alternative  approach  to  the  reconstruction  is  to  apply  the  stencil  shown  in  Equations  12  and  13  to 
the  characteristic  variables  Win  the  computational  space,  where  W  =  R~'Q .  Once  W^+1/2  and 

Wm+i/2 are  obtained,  the  independent  variables  Q^+1/2  and  Qj^+1/2  can  be  recovered  by  pre- 
multipling  W*+1/2  and  W^+]/2  with  the  right  eignvector  R  computed  at  the  interface  m+1/2.  This 

reconstruction  process,  identified  here  as  weighted  MUSCL  here,  requires  higher  operation  count, 
but  it  takes  the  coordinate  non-uniformity  and  stretching  into  account.  The  scheme  described 
here,  has  been  used  to  compute  the  flow  in  linear  cascade  with  tip  clearance  [12].  Similarly,  a  fifth 
order  stencil,  suggested  by  Rai  [13]  and  used  by  Rogers  and  Kwak  [14],  is  derived  and  applied  to 
the  characteristic  variables  here.  That  is 

Wm\1/2  =  Wm+1  +  ^  [+2Am+5/2  W  - 1 1  Ara+3/2W  -  24Am+1/2  W  +  3Am_1/2W]  (14) 

and 

Wm\1/2  =  Wm  +^[-2Am.3/2W  +  llAm„1/2W  +  24Am+1/2W-3Am+3/2W]  (15) 

The  weighted  MUSCL  schemes  described  here  are  implemented  in  UNCLE  code  for  present 
calculations. 

The  numerical  viscous  fluxes  defined  in  Equation  1  can  be  obtained  through  central 
differencing  technique.  The  differencing  equations  that  represent  the  Navier-Stokes  equation 
formulated  in  conservation  form  can  then  be  assembled.  The  system  of  equations  is  solved  with 
Newton  iteration  method.  Standard  ks  type  turbulence  equation  is  used  for  turbulence  simulation. 

GRID  GENERATION 

The  inflow  to  the  rotor  is  assumed  axisymmetrical;  therefore,  only  one  blade  passage, 
bounded  between  two  adjacent  blades,  with  appropriated  boundary  condition  is  sufficient  for  the 
numerical  simulation.  The  water  tunnel,  as  shown  in  Figure  1,  is  an  open  jet  type,  with  jet 
discharging  into  a  plenum  from  the  tunnel  nozzle  downstream  of  the  rotor.  At  current  simulation, 
the  plenum  is  not  modeled,  and  the  tunnel  is  assumed  to  extend  downstream  with  constant  radius, 
with  no-slip  boundary  condition  applied  on  its  wall.  In  the  blade  passage,  a  half  C-type  grid  is 
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used  to  wrap  around  the  pressure  side  and  suction  side  of  the  blade,  covering  the  volume 
immediately  around  the  blade  surface.  The  rest  of  the  volume,  extended  4.5  rotor  radius  both 
upstream  and  downstream  from  the  stacking  line,  including  the  tip  gap,  is  represented  by  the  H- 
type  grid.  The  C-type  grid  wrapping  around  the  blade  is  to  provide  a  better  grid  quality  near  the 
leading  edge.  The  first  grid  point  off  the  blade  surface  is  about  five  millionth  of  the  blade  radius, 
at  the  testing  Reynolds  number  (3.1  millions,  based  on  the  tunnel  speed  and  blade  radius),  the 
corresponding  y+  is  about  2.  Based  on  the  topology  just  described,  two  grids  with  different  nodal 
numbers  are  generated.  The  coarse  grid  has  about  2.8  million  nodes  and  the  fine  grid  has  about 
5.5  million  nodes.  The  grid  system  is  then  divided  into  ten  connecting  blocks  for  multi-processors 
processing. 


COMPUTATIONAL  RESULTS 

The  computations  were  performed  on  SGI  3400  multi-processing  machines.  The 
solutions  were  interpolated  onto  the  measurement  planes  for  the  purpose  of  comparisons  with  the 
experimental  data.  The  velocity  profile  shown  in  Figure  2  was  used  as  upstream  inflow  for 
present  computation.  The  reference  static  pressure  was  fixed  to  zero  at  far  field  downstream. 

Figure  4  shows  the  formation  of  the  tip  vortices  from  the  computation  based  on  the 
coarse  grid.  The  flow  enters  the  tip  gap  region  from  the  pressure  side  of  the  blade  and  vorticity  is 
generated  and  convected  in  the  gap  region.  As  the  flow  exits  to  the  pressure  side,  the  roll-up 
begins.  As  the  flow  exits  near  the  leading  edge,  the  vortex  rolls  up  and  forms  the  core  of  the 
leakage  vortex.  The  flow  exits  around  the  mid-chord  and  after  rolls  up  and  wraps  around  the  core. 
At  the  trailing  edge,  the  flow  from  the  tip  gap  merges  with  the  flow  from  the  pressure  side  and  the 
flow  from  the  suction  side  and  forms  the  trailing  edge  vortex.  The  trailing  edge  vortex  then  wraps 
around  the  leakage  vortex.  The  region  of  the  depression  of  the  static  pressure  coincides  with  the 
passage  of  the  vortex  core. 

Figure  5  shows  the  contours  of  the  computed  static  pressure  at  various  x-r  planes  behind 
the  trailing  edge.  There  are  three  sets  of  images  in  Figure  5.  The  image  (a)  is  based  on  the  coarse 
grid  with  standard  MUSCL  scheme;  the  images  (b)  and  (c)  are  based  on  fine  grid  with  third  and 
fifth  order  weighted  MUSCL  schemes,  respectively.  Comparison  of  the  three  images  indicates 
that  the  pressure  distribution  is  less  diffusive  when  finer  grid  or  higher  order  differencing  scheme 
is  employed. 

One  way  to  trace  the  vortex  in  the  three-dimensional  flow  field  is  through  the  display  of 
helicity,  which  is  the  projection  of  the  velocity  vector  V  onto  the  vorticity  vector  m  ,  that 
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is  H  =  V  •  G7 .  Figure  6  shows  the  helicity  downstream  of  tip  trailing  edge.  The  image  (a)  is  the 
computed  helicity  with  the  third  order  weighted  MUSCL  scheme  based  on  the  fine  grid,  the 
image  (b)  is  the  helicity  derived  from  the  measurement.  Based  on  the  intensity  of  the  helicity,  it 
can  be  seen  from  both  the  computation  and  measurement  that  at  a  distance  about  10%  chord 
length  downstream  of  the  trailing  edge,  the  trailing  edge  vortex  crosses  over  and  wraps  around  the 
leakage  vortex.  At  that  location,  the  computed  helicity  is  slightly  less  than  measured.  The 
contours  displayed  in  Figure  6  suggest  that  the  computed  trajectory  of  the  vortex  center  is  slightly 
further  away  from  the  suction  side  than  measured. 

Figure  7  to  Figure  12  show  the  tangential  velocity  vectors  and  streamwise  vorticity 
contours  on  the  x-r  planes  at  a  distance  of  -0.0108,  0.0596,  0.0920,  0.1245  and  0.3032  chord 
length  from  the  trailing  edge,  respectively.  The  tangential  velocity  is  obtained  by  projecting  the 
velocity  vector  onto  the  x-r  plane  and  the  streamwise  vorticity  is  obtained  by  dividing  the  helicity 
with  the  magnitude  of  the  velocity.  In  each  of  those  figures,  two  images  are  presented.  The  image 
on  the  left  is  derived  from  the  experimental  data  and  the  image  on  the  right  is  from  the 
computational  results.  Figure  7  shows  that  at  the  immediately  upstream  of  the  trailing  edge,  the 
vorticity  is  concentrated  in  the  leakage  vortex  and  in  the  shear  layer  that  exits  the  gap.  The 
vorticity  in  the  shear  layer  feeds  into  the  leakage  vortex  continuously;  in  the  meantime,  the 
viscous  diffusion  takes  place.  Those  two  competing  elements  determine  the  strength  of  the 
vorticity  in  the  vortex.  The  locations  of  the  computed  and  the  measured  centers  of  vortex, 
estimated  by  the  orientation  of  the  tangential  velocity  vectors  and  the  intensity  of  the  vorticity, 
match  quite  well.  However,  the  intensity  of  the  computed  streamwise  vorticity  is  lower  than  the 
measured.  The  discrepancy  may  be  caused  by  the  grid  resolution  in  the  region  of  vortex  passage. 
In  a  earlier  vortex  simulation,  it  was  found  that  to  resolve  the  cross  section  of  the  vortex  core 
accurately,  a  20  by  20  mesh  is  needed.  At  present,  the  average  mesh  size  in  core  region  is  about 
10  by  10. 

Figure  8  shows  that  at  a  distance  of  6%  chord  length  downstream  of  the  trailing  edge,  in 
addition  to  the  leakage  vortex,  a  trailing  edge  vortex  appears.  Its  present  can  be  seen  by  the 
intense  vorticity  concentration  above  the  leakage  vortex.  Figure  9  and  Figure  10  show  that  the 
trailing  edge  vortex  wraps  around  the  leakage  vortex,  as  it  convects  further  downstream.  At  this 
stage,  the  predicted  trajectory  of  the  vortex  center  becomes  further  away  from  the  blade  suction 
surface  than  the  one  measured.  Figure  1 1  and  Figure  12  show  that  at  further  downstream  from 
trailing  edge,  at  s/C  =  0.3010  and  s/C  =  0.4037,  the  trailing  edge  vortex  becomes  very  weak  and 
undistinguishable. 

Figure  13  to  Figure  17  show  the  computed  tangential  velocity  vectors  and  static  pressure 
contours  on  the  x-r  planes  at  a  distance  of  0.0%,  10%,  20%,  30%  and  40%  chord  length 
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downstream  from  the  trailing  edge.  The  results  are  based  on:  (a)  coarse  grid  with  the  standard 
MUSCL  scheme,  (b)  fine  grid  with  the  standard  MUSCL  scheme,  (c)  fine  grid  with  the  weighted 
MUSCL  scheme,  and  (d)  fine  grid  with  the  fifth-order  MUSCL  scheme.  Through  the  pressure 
contours,  Figure  14  indicates  that  at  s/C=0.1026  the  trailing  edge  vortex  is  better  resolved  with 
weighted  MUSCL  scheme,  especially  with  the  fifth-order  accuracy. 

Figure  18  shows  the  computed  static  pressure  distribution  in  the  vortex  core.  Results 
presented  are  based  on  both  coarse  and  fine  grids,  with  MUSCL  and  weighted  MUSCL  schemes. 
For  the  purpose  of  comparison,  Iowa-Kim’s  result  is  also  presented  [5].  The  reference  pressure 
for  present  calculations  is  set  at  zero  at  downstream  exit.  This  is  compatible  with  Kim’s  approach 
[5].  For  all  the  calculations,  the  minimum  static  pressure  within  the  flow  field  occurs  at  a  distance 
about  1 1%  chord  length  downstream  of  the  trailing  edge  where  the  trailing  edge  vortex  is 
crossing  over  the  leakage  vortex  and  begins  the  wrapping  around  process.  Overall,  at  a  given 
location,  the  fine  grid  calculation  produces  lower  vortex  core  pressure  than  the  coarse  grid 
calculation;  the  weighted  MUSCL  scheme  produces  lower  vortex  core  pressure  than  the  standard 
MUSCL  scheme,  and  the  fifth  order  differencing  scheme  produces  even  lower  level  of  pressure. 
The  computed  minimum  static  pressure  on  blade  surface  is  about  -4.0,  and  is  located  on  the 
suction  side  of  the  blade  tip,  about  80%  chord  length  from  the  leading  edge.  This  concurs  with 
Kim’s  finding  [5], 

In  general,  at  a  given  x-r  plan,  the  location  of  minimum  pressure  coincides  with  that  of 
the  vortex  center.  With  the  same  differencing  scheme,  the  solution  of  the  fine  grid  gives  lower 
vortex  core  pressure  than  the  solution  of  the  coarse  grid.  With  the  same  grid,  the  solution  of  the 
weighted  scheme  gives  lower  vortex  core  pressure  than  the  solution  of  the  unweighted  scheme. 
The  fifth-order  differencing  scheme  gives  the  lowest  vortex  core  pressure.  The  pressure 
differentiations  among  the  computations  mentioned  above  take  place  only  within  the  vortex  core 
where  the  gradients  of  velocity  and  pressure  are  sharp.  In  the  region  away  from  the  vortex  core, 
the  effects  of  the  grid  density  and  the  order  of  the  differencing  schemes  are  less  significant. 

DISCUSSION  AND  CONCLUSION 

A  successful  numerical  simulation  of  vortex  flow  requires  higher  order  differencing 
scheme  and  dense  grid  support.  At  present,  to  simulate  the  vortex  flow  about  the  ducted-rotor 
P5206,  two  sets  of  grid  with  different  nodal  numbers  (2.8  million  and  5.5  million  nodes)  are 
utilized  for  grid  sensitivity  study,  and  a  fifth-order  weighted  MUSCL  scheme  is  derived  and 
implemented  to  reduce  the  numerical  diffusion  in  the  vortex  core.  Numerical  solutions  show  that 
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higher  order  differencing  scheme  indeed  gives  sharper  and  less  diffusive  static  pressure  and 
streamwise  vorticity  distributions  within  the  vortex  core. 

Results  from  the  fine  grid  and  fifth-order  solution  show  that  the  minimum  pressure  in  the 
flow  field  occurs  in  the  vortex  core  at  about  11%  chord  length  downstream  of  the  tip  trailing 
edge.  This  is  the  location  where  the  trailing  edge  vortex  crossing  over  the  leakage  vortex.  At  that 
location,  the  pressure  coefficient  Cpis  about  -10.22.  According  to  the  observation  of  cavitation 
inception,  at  the  rotor  operating  condition,  the  cavitation  index  C;  is  about  11.  The  location  of 
inception  is  about  half  a  chord  length  downstream  of  the  trailing  edge.  At  that  location,  the 
measured  helicity  or  streamwise  vorticity  do  not  show  any  unusually  high  values  to  indicate  the 
possibility  of  the  existence  of  low  static  pressure.  The  cause  of  the  discrepancy  may  attribute  to 
the  unsteadiness  of  the  vertical  flow.  The  unsteady  nature  of  the  trailing  vortices  has  been 
reported  by  Green  [15,16,17],  and  Green  and  Acosta  [18].  It  has  been  observed  that  the  flow  in 
the  vortex  core  is  highly  unsteady  and  the  core  pressure  fluctuations  are  a  remarkable  one  to  two 
times  the  freestream  dynamic  pressure.  Outside  of  the  core,  the  flow  unsteadiness  diminishes 
rapidly  and  the  measured  flow  unsteadiness  is  largely  due  to  the  freestream  turbulence.  It  has  also 
been  observed  that  the  region  that  was  0.3  to  2  chord  length  downstream  from  the  trailing  edge, 
was  the  most  unsteady  portion  of  the  vortex  core.  Due  to  the  large  static  pressure  fluctuation,  the 
cavitation  inception  occurred  at  mean  core  pressure  above  the  vapor  pressure. 

In  the  investigation  of  turbulent  shear  cavitation,  O’Hem  [19]  found  that  in  the  vortex 
core,  the  pressure  fluctuations  had  positive  and  negative  peaks  as  high  as  three  times  the 
freestream  dynamic  pressure,  sufficient  to  explain  cavitation  inception  at  high  values  of  the 
inception  index.  In  summary,  the  flow  in  the  vortex  core  is  highly  unsteady,  especially  in  the 
region  0.3  to  2  times  chord  length  downstream  of  the  trailing  edge.  It  is  highly  likely  that  the 
inception  may  occur  at  mean  pressure  higher  than  the  vapor  pressure  in  that  region.  At  present 
computations,  the  minimum  static  pressure  occurs  at  a  distance  about  0.1 1  chord  length 
downstream  of  the  trail  edge,  where  the  flow  is  mostly  steady  and  Reynolds  number  independent. 
At  a  distance  0.5  chord  length  downstream  of  the  trailing  edge,  the  flow  is  mostly  unsteady,  and 
the  predicted  mean  static  pressure  is  -7.0.  A  negative  pressure  fluctuation  of  -4.0  is  sufficient  to 
cause  the  inception.  Such  a  magnitude  of  fluctuation  is  likely  to  occur  according  to  Green’s 
studies[15,16,17,18]. 
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Figure  7:  Computed  and  measured  tangential  velocity  vectors  and  streamwise  vorticity 

contours  at  s/C=-0.0108. 
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Figure  8:  Computed  and  measured  tangential  velocity  vectors  and  streamwise  vorticity 

contours  at  s/C  =  0.0596. 
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Figure  9:  Computed  and  measured  tangential  velocity  vectors  and  streamwise  vorticity 

contours  at  s/C=0.0920. 
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Figure  10:  Computed  and  measured  tangential  velocity  vectors  and  streamwise  vorticity 

contours  at  s/C=0.1245. 
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Figure  11:  Computed  and  measured  tangential  velocity  vectors  and  streamwise  vorticity 

contours  at  s/C=0.3032. 
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Figure  12:  Computed  and  measured  tangential  velocity  vectors  and  steamwise  vorticity 

contours  at  s/C=0.4006. 


NSWCCD-50-TR-2004/041 


19 


NSWCCD-50-TR-2004/041 


(a)  Coarse  grid  with  MUSCL  (b)  Fine  grid  with  MUSCL 


* 


(c)  Fine  grid  with  weighted  MUSCL 


Figure  14:  Computed  tangential  velocity  vectors  and  static  pressure  contours  at 

s/C=0.1026. 
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(a)  Coarse  grid  with  MUSCL 


(b)  Fine  grid  with  MUSCL 


(c)  Fine  grid  with  weighted  MUSCL  (d)  Fine  grid  with  5th-order  weighted  MUSCL 


Figure  16:  Computed  tangential  velocity  vectors  and  static  pressure  contours  at 

s/C=0.3031. 
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